subroutine allocate_camera_arrays()

#ifdef USE_CAM_FILE

  use camera_module
  implicit none
 
  if(allocated(cam_Frames)) deallocate(cam_Frames)
  allocate(cam_Frames(NumFramesCamFile))
  if(.NOT.(allocated(cam_Frames))) call abort("Failed to allocate cam_Frames")

  if(allocated(cam_x)) deallocate(cam_x)
  allocate(cam_x(NumFramesCamFile))
  if(.NOT.(allocated(cam_x))) call abort("Failed to allocate cam_x")

  if(allocated(cam_y)) deallocate(cam_y)
  allocate(cam_y(NumFramesCamFile))
  if(.NOT.(allocated(cam_y))) call abort("Failed to allocate cam_y")

  if(allocated(cam_z)) deallocate(cam_z)
  allocate(cam_z(NumFramesCamFile))
  if(.NOT.(allocated(cam_z))) call abort("Failed to allocate cam_z")

#ifdef VELOCITY
  if(allocated(cam_vx)) deallocate(cam_vx)
  allocate(cam_vx(NumFramesCamFile))
  if(.NOT.(allocated(cam_vx))) call abort("Failed to allocate cam_vx")

  if(allocated(cam_vy)) deallocate(cam_vy)
  allocate(cam_vy(NumFramesCamFile))
  if(.NOT.(allocated(cam_vy))) call abort("Failed to allocate cam_vy")

  if(allocated(cam_vz)) deallocate(cam_vz)
  allocate(cam_vz(NumFramesCamFile))
  if(.NOT.(allocated(cam_vz))) call abort("Failed to allocate cam_vz")
#endif

  if(allocated(cam_Lx)) deallocate(cam_Lx)
  allocate(cam_Lx(NumFramesCamFile))
  if(.NOT.(allocated(cam_Lx))) call abort("Failed to allocate cam_Lx")

  if(allocated(cam_Ly)) deallocate(cam_Ly)
  allocate(cam_Ly(NumFramesCamFile))
  if(.NOT.(allocated(cam_Ly))) call abort("Failed to allocate cam_Ly")

  if(allocated(cam_Lz)) deallocate(cam_Lz)
  allocate(cam_Lz(NumFramesCamFile))
  if(.NOT.(allocated(cam_Lz))) call abort("Failed to allocate cam_Lz")

#ifdef USE_EIGENVECTORS
  if(allocated(cam_xprime)) deallocate(cam_xprime)
  allocate(cam_xprime(3,NumFramesCamFile))
  if(.NOT.(allocated(cam_xprime))) call abort("Failed to allocate cam_xprime")

  if(allocated(cam_yprime)) deallocate(cam_yprime)
  allocate(cam_yprime(3,NumFramesCamFile))
  if(.NOT.(allocated(cam_yprime))) call abort("Failed to allocate cam_yprime")

  if(allocated(cam_zprime)) deallocate(cam_zprime)
  allocate(cam_zprime(3,NumFramesCamFile))
  if(.NOT.(allocated(cam_zprime))) call abort("Failed to allocate cam_zprime")
#else
  if(allocated(cam_theta)) deallocate(cam_theta)
  allocate(cam_theta(NumFramesCamFile))
  if(.NOT.(allocated(cam_theta))) call abort("Failed to allocate cam_theta")

  if(allocated(cam_phi)) deallocate(cam_phi)
  allocate(cam_phi(NumFramesCamFile))
  if(.NOT.(allocated(cam_phi))) call abort("Failed to allocate cam_phi")

  if(allocated(cam_psi)) deallocate(cam_psi)
  allocate(cam_psi(NumFramesCamFile))
  if(.NOT.(allocated(cam_psi))) call abort("Failed to allocate cam_psi")
#endif
#endif

end subroutine allocate_camera_arrays

! --------------------------------------------------

subroutine allocate_index()
  use my_mpi
  use particles, only : indx, NumPart_ByPE
  implicit none
  
  ! Allocates memory for the sort index array of the
  ! local particles
  if(allocated(indx)) deallocate(indx)
  allocate(indx(NumPart_ByPE(MyPE)))
  if(.not.(allocated(indx))) call abort(' -failed to allocate indx')

end subroutine allocate_index

! --------------------------------------------------

subroutine deallocate_index()
  use particles, only : indx
  implicit none
  
  if(allocated(indx)) deallocate(indx)

end subroutine deallocate_index

! --------------------------------------------------

subroutine allocate_exchange_counters()
  use my_mpi, only : NumPEs
  use exchange
  implicit none

  ! Allocates memory for the import and export arrays,  an
  ! offset array for each, and the ExchangeMatrix
  if(allocated(NumImportPE)) deallocate(NumImportPE)
  allocate(NumImportPE(NumPEs))
  if(.not.(allocated(NumImportPE))) call abort(' -failed to allocate NumImportPE')

  if(allocated(NumExportPE)) deallocate(NumExportPE)
  allocate(NumExportPE(NumPEs))
  if(.not.(allocated(NumExportPE))) call abort(' -failed to allocate NumExportPE')

  if(allocated(ImportOffset)) deallocate(ImportOffset)
  allocate(ImportOffset(NumPEs))
  if(.not.(allocated(ImportOffset))) call abort(' -failed to allocate ImportOffset')

  if(allocated(ExportOffset)) deallocate(ExportOffset)
  allocate(ExportOffset(NumPEs))
  if(.not.(allocated(ExportOffset))) call abort(' -failed to allocate ExportOffset')

  if(allocated(ExchangeMatrix)) deallocate(ExchangeMatrix)
  allocate(ExchangeMatrix(NumPEs*NumPEs))
  if(.not.(allocated(ExchangeMatrix))) call abort(' -failed to allocate ExchangeMatrix')

end subroutine allocate_exchange_counters

! --------------------------------------------------

subroutine deallocate_exchange_counters()
  use exchange
  implicit none

  if(allocated(NumImportPE)) deallocate(NumImportPE)
  if(allocated(NumExportPE)) deallocate(NumExportPE)
  if(allocated(ImportOffset)) deallocate(ImportOffset)
  if(allocated(ExportOffset)) deallocate(ExportOffset)
  if(allocated(ExchangeMatrix)) deallocate(ExchangeMatrix)

end subroutine deallocate_exchange_counters

! --------------------------------------------------

subroutine allocate_destination_array
  use my_mpi
  use particles, only : NumPart_ByPE
  use exchange
  implicit none

  if(allocated(DestinationPE)) deallocate(DestinationPE)
  allocate(DestinationPE(NumPart_ByPE(MyPE)))
  if(.not.(allocated(DestinationPE))) call abort(' -failed to allocate DestinationPE')

end subroutine allocate_destination_array

! --------------------------------------------------

subroutine deallocate_destination_array
  use exchange
  implicit none

  if(allocated(DestinationPE)) deallocate(DestinationPE)

end subroutine deallocate_destination_array

! --------------------------------------------------

subroutine allocate_import_buffer(itype,nbytes)
  use logs
  use my_mpi
  use exchange, only : ImportBuf, NumImport, ImportBuf_Int
  implicit none
  integer*4 :: itype, nbytes
  integer*4 :: i,ierr,position
  integer*8 :: nbytes_tot

  ! Ensure at least a single byte is allocated, otherwise
  ! code fails during mpi_sendrecv
  nbytes_tot = dble(nbytes)*dble(NumImport)

  if(nbytes_tot < 0) then 
     write(logunit,*)' in allocate_import_buffer: nbytes_tot = ',nbytes_tot
     write(logunit,*)' nbyte per particle = ',nbytes
     write(logunit,*)' number of particles = ',NumImport
     call abort("Overflow in import buffer")
  endif

  if(allocated(ImportBuf)) deallocate(ImportBuf)
  allocate(ImportBuf(max(nbytes_tot,1)))
  if(.not.(allocated(ImportBuf))) call abort(' -failed to allocate ImportBuf')

end subroutine allocate_import_buffer

! --------------------------------------------------

subroutine deallocate_import_buffer()
  use exchange, only : ImportBuf
  implicit none

  if(allocated(ImportBuf)) deallocate(ImportBuf)  
  
end subroutine deallocate_import_buffer

! --------------------------------------------------

subroutine allocate_export_buffer(itype,nbytes)
  use logs
  use my_mpi, only : NumPEs
  use exchange, only : ExportBuf, NumExport
  implicit none
  integer*4 :: itype, nbytes
  integer*8 :: nbytes_tot

  ! nbytes is per particle, nbytes_tot is total needed for the buffer
  ! We seem to need an extra byte when using mpi_pack since the
  ! position indicator rolls on by 1 after final variable is packed
  nbytes_tot = (nbytes*NumExport)  +1

  if(nbytes_tot < 0) then 
     write(logunit,*)' in allocate_export_buffer: nbytes_tot = ',nbytes_tot
     call abort("Overflow in export buffer")
  endif

  if(allocated(ExportBuf)) deallocate(ExportBuf)
  allocate(ExportBuf(nbytes_tot))
  if(.not.(allocated(ExportBuf))) call abort(' -failed to allocate ExportBuf')

end subroutine allocate_export_buffer

! --------------------------------------------------

subroutine deallocate_export_buffer()
  use exchange, only : ExportBuf
  implicit none

  if(allocated(ExportBuf)) deallocate(ExportBuf)
  
end subroutine deallocate_export_buffer

! --------------------------------------------------

subroutine allocate_local_grids(itype)
  use logs
  use parameter_module
  use grid
  use smoothing_variables, only : nvar_by_type
  implicit none
  integer*4 :: itype
  
  ! Count many grids we need to allocate
  ! If using SPH we need an extra grid for the normalisation
  ngrids = nvar_by_type(itype)
    
#ifdef CUBE
  if(allocated(ValueQuantity)) deallocate(ValueQuantity)
  allocate(ValueQuantity(ngrids,params%nx,params%ny,nzlocal)) 
  ValueQuantity(:,:,:,:) = 0.0e0
#else
  if(allocated(ValueQuantity)) deallocate(ValueQuantity)
  allocate(ValueQuantity(ngrids,params%nx,nzlocal)) 
  ValueQuantity(:,:,:) = 0.0e0
#endif	

  if(allocated(ValueQuantity))then
     write(mystring,'(F7.1)')real(4.0*size(ValueQuantity))/(1024.0*1024.0)	     
     write(logunit,*)'   -allocated '//trim(adjustl(mystring))//'Mb for ValueQuantity grid'
  else
     call abort("failed to allocate ValueQuantity grid")
  endif

  ! Also allocate Value grid.
#ifdef CUBE
  if(allocated(Value)) deallocate(Value)
  allocate(Value(params%nx,params%ny,nzlocal)) 
  Value(:,:,:) = 0.0e0
#else
  if(allocated(Value)) deallocate(Value)
  allocate(Value(params%nx,nzlocal)) 
  Value(:,:) = 0.0e0
#endif

  if(allocated(Value))then
     write(mystring,'(F7.1)')real(4.0*size(Value))/(1024.0*1024.0)	     
     write(logunit,*)'   -allocated '//trim(adjustl(mystring))//'Mb for Value grid'
  else
     call abort("failed to allocate Value grid")
  endif

  ! Also do the min,max,mean arrays
  If(allocated(ValueQuantityMin)) deallocate(ValueQuantityMin)
  allocate(ValueQuantityMin(ngrids))  
  if(allocated(ValueQuantityMax)) deallocate(ValueQuantityMax)
  allocate(ValueQuantityMax(ngrids))  
  if(allocated(ValueQuantitySum)) deallocate(ValueQuantitySum)
  allocate(ValueQuantitySum(ngrids))  

end subroutine allocate_local_grids

! --------------------------------------------------

subroutine deallocate_local_grids
  use grid, only : Value, ValueQuantity
  implicit none

  if(allocated(Value))         deallocate(Value)
  if(allocated(ValueQuantity)) deallocate(ValueQuantity)

end subroutine deallocate_local_grids

! --------------------------------------------------

subroutine allocate_read_particles(itype)
  use my_mpi
  use logs
  use particles
  use smoothing_variables
  implicit  none
  integer*4 :: nb
  integer*4 :: itype

  if(allocated(readCoordinates)) deallocate(readCoordinates)
  allocate(readCoordinates(3,NumPart_ByPE(MyPE)))
  if(.not.(allocated(readCoordinates))) call abort("Failed to allocate readCoordinates")
  nb = 3*NumPart_ByPE(MyPE)

#ifdef VELOCITY
     if(allocated(readVelocity)) deallocate(readVelocity)
     allocate(readVelocity(3,NumPart_ByPE(MyPE)))
     if(.not.(allocated(readVelocity))) call abort("Failed to allocate readVelocity")
     nb = 3*NumPart_ByPE(MyPE)
#endif

  if(allocated(readMass)) deallocate(readMass)
  allocate(readMass(NumPart_ByPE(MyPE)))
  if(.not.(allocated(readMass))) call abort("Failed to allocate readMass")
  nb = nb + NumPart_ByPE(MyPE)

  if(allocated(readSmoothingLength)) deallocate(readSmoothingLength)
  allocate(readSmoothingLength(NumPart_ByPE(MyPE)))
  if(.not.(allocated(readSmoothingLength))) call abort("Failed to allocate readSmoothingLength")
  nb = nb + NumPart_ByPE(MyPE)

#ifndef FORCE_EXTERNAL_HSML
  if(itype == 0) then
#endif
     if(allocated(readParticleIDs)) deallocate(readParticleIDs)
     allocate(readParticleIDs(NumPart_ByPE(MyPE)))
     if(.not.(allocated(readParticleIDs))) call abort("Failed to allocate readParticleIDs")
     nb = nb + NumPart_ByPE(MyPE)
#ifndef FORCE_EXTERNAL_HSML
endif
#endif

#ifdef FORCE_EOS
  if(itype == 0) then
     if(allocated(readEoS)) deallocate(readEoS)
     allocate(readEoS(NumPart_ByPE(MyPE)))
     if(.not.(allocated(readEoS))) call abort("Failed to allocate readEoS")
     nb = nb + NumPart_ByPE(MyPE)
  endif
#endif

  if(allocated(readSmoothVar)) deallocate(readSmoothVar)
  allocate(readSmoothVar(nvar_by_type(itype),NumPart_ByPE(MyPE)))
  if(.not.(allocated(readSmoothVar))) call abort("Failed to allocate readSmoothVar")
  nb = nb + nvar_by_type(itype)*NumPart_ByPE(MyPE)

  write(mystring,'(F7.1)')4.0e0*real(nb)/(1024.0*1024.0)
  write(logunit,*)'   -allocated '//trim(adjustl(mystring))//'Mb for i/o phase particle storage'
  call flush(logunit)

end subroutine allocate_read_particles

! -------------------------------------------------

subroutine deallocate_read_particles
  use particles
  implicit  none

  if(allocated(readCoordinates)) deallocate(readCoordinates)
#ifdef VELOCITY
  if(allocated(readVelocity)) deallocate(readVelocity)
#endif
  if(allocated(readMass)) deallocate(readMass)  
  if(allocated(readSmoothingLength)) deallocate(readSmoothingLength)
  if(allocated(readParticleIDs)) deallocate(readParticleIDs)
  if(allocated(readEoS)) deallocate(readEoS)
  if(allocated(readSmoothVar)) deallocate(readSmoothVar)
  
end subroutine deallocate_read_particles

! --------------------------------------------------

subroutine allocate_particles(itype)
  use my_mpi
  use logs
  use particles
  use smoothing_variables, only : nvar_by_type
  implicit  none
  integer*4 :: nb
  integer*4 :: itype

  if(allocated(Coordinates)) deallocate(Coordinates)
  allocate(Coordinates(3,NumPart_ByPE(MyPE))) ; Coordinates(:,:) = 0.0e0
  if(.not.(allocated(Coordinates))) call abort("Failed to allocate Coordinates")
  nb = 3*NumPart_ByPE(MyPE)

  if(allocated(Mass)) deallocate(Mass)
  allocate(Mass(NumPart_ByPE(MyPE))) ; Mass(:) = 0.0e0
  if(.not.(allocated(Mass))) call abort("Failed to allocate Mass")
  nb = nb + NumPart_ByPE(MyPE)

  if(allocated(SmoothingLength)) deallocate(SmoothingLength)
  allocate(SmoothingLength(NumPart_ByPE(MyPE))) ; SmoothingLength(:) = 0.0e0
  if(.not.(allocated(SmoothingLength))) call abort("Failed to allocate SmoothingLength")
  nb = nb + NumPart_ByPE(MyPE)

#ifdef FORCE_EOS
  if(itype == 0) then
     if(allocated(EoS)) deallocate(EoS)
     allocate(EoS(NumPart_ByPE(MyPE))) ; EoS(:) = 0.0e0
     if(.not.(allocated(EoS))) call abort("Failed to allocate EoS")
     nb = nb + NumPart_ByPE(MyPE)
  endif
#endif

  if(allocated(SmoothVar)) deallocate(SmoothVar)
  allocate(SmoothVar(nvar_by_type(itype),NumPart_ByPE(MyPE))) ; SmoothVar(:,:) = 0.0e0
  if(.not.(allocated(SmoothVar))) call abort("Failed to allocate SmoothVar")
  nb = nb + nvar_by_type(itype)*NumPart_ByPE(MyPE)

  write(mystring,'(F7.1)')4.0e0*real(nb)/(1024.0*1024.0)
  write(logunit,*)'   -allocated '//trim(adjustl(mystring))//'Mb for particle storage'
  call flush(logunit)

end subroutine allocate_particles

! --------------------------------------------------

subroutine allocate_ghost_particles(itype)
  use my_mpi
  use logs
  use smoothing_variables, only : nvar_by_type
  use particles
  implicit  none
  integer*4 :: nb
  integer*4 :: itype

  if(allocated(Ghost_Coordinates)) deallocate(Ghost_Coordinates)
  allocate(Ghost_Coordinates(3,NumGhost_ByPE(MyPE))) ; Ghost_Coordinates(:,:) = 0.0e0
  if(.not.(allocated(Ghost_Coordinates))) call abort("Failed to allocate Ghost_Coordinates")
  nb = 3*NumGhost_ByPE(MyPE)

  if(allocated(Ghost_Mass)) deallocate(Ghost_Mass)
  allocate(Ghost_Mass(NumGhost_ByPE(MyPE))) ; Ghost_Mass(:) = 0.0e0
  if(.not.(allocated(Ghost_Mass))) call abort("Failed to allocate Ghost_Mass")
  nb = nb + NumGhost_ByPE(MyPE)

  if(allocated(Ghost_SmoothingLength)) deallocate(Ghost_SmoothingLength)
  allocate(Ghost_SmoothingLength(NumGhost_ByPE(MyPE))) ; Ghost_SmoothingLength(:) = 0.0e0
  if(.not.(allocated(Ghost_SmoothingLength))) call abort("Failed to allocate Ghost_SmoothingLength")
  nb = nb + NumGhost_ByPE(MyPE)

#ifdef FORCE_EOS
  if(itype == 0) then
     if(allocated(Ghost_EoS)) deallocate(Ghost_EoS)
     allocate(Ghost_EoS(NumGhost_ByPE(MyPE))) ; Ghost_EoS(:) = 0.0e0
     if(.not.(allocated(Ghost_EoS))) call abort("Failed to allocate Ghost_EoS")
     nb = nb + NumPart_ByPE(MyPE)
  endif
#endif

  if(allocated(Ghost_SmoothVar)) deallocate(Ghost_SmoothVar)
  allocate(Ghost_SmoothVar(nvar_by_type(itype),NumGhost_ByPE(MyPE))) ; Ghost_SmoothVar(:,:) = 0.0e0
  if(.not.(allocated(Ghost_SmoothVar))) call abort("Failed to allocate Ghost_SmoothVar")
  nb = nb + NumGhost_ByPE(MyPE)
  
  write(mystring,'(F7.1)')4.0e0*real(nb)/(1024.0*1024.0)
  write(logunit,*)'   -allocated '//trim(adjustl(mystring))//'Mb for ghost particle storage'
  call flush(logunit)

end subroutine allocate_ghost_particles

! ------------------------------------------------

subroutine deallocate_ghost_particles
  use particles
  implicit none

  if(allocated(Ghost_Coordinates)) deallocate(Ghost_Coordinates)
  if(allocated(Ghost_Mass)) deallocate(Ghost_Mass)
  if(allocated(Ghost_SmoothingLength)) deallocate(Ghost_SmoothingLength)
  if(allocated(Ghost_SmoothVar)) deallocate(Ghost_SmoothVar)

end subroutine deallocate_ghost_particles

! ------------------------------------------------

subroutine allocate_hsmlfile()
  use logs
  use header
  use hsmlfile
  implicit none

  if(allocated(hsml%ParticleIDs)) deallocate(hsml%ParticleIDs)
  allocate(hsml%ParticleIDs(NumPart_Total(1))) ; hsml%ParticleIDs(:) = 0
  if(.not.(allocated(hsml%ParticleIDs))) call abort("Failed to allocate hsml%ParticleIDs")

  if(allocated(hsml%SmoothingLength)) deallocate(hsml%SmoothingLength)
  allocate(hsml%SmoothingLength(NumPart_Total(1))) ; hsml%SmoothingLength(:) = 0.0e0
  if(.not.(allocated(hsml%SmoothingLength))) call abort("Failed to allocate hsml%SmoothingLength")
  

end subroutine allocate_hsmlfile
